home *** CD-ROM | disk | FTP | other *** search
- %
- % "simplex.t" solves a linear programming problem
- % using the simplex method
- %
- % Sample program for the T Interpreter by:
- %
- % Stephen R. Schmitt
- % 962 Depot Road
- % Boxborough, MA 01719
- %
-
- const N : int := 3 % number of variables
- const M : int := 4 % number of conditions
-
- var Sx, Tx : array[M+1,N+1] of real
- var basis : array[M+1] of int
- var non_basis : array[N+1] of int
-
- program
-
- var row, col : int
- label program_exit :
-
- initialize
-
- loop
-
- show_array
-
- col := pivot_col
- exit when col = 0 % all done
-
- row := pivot_row( col )
- if row = 0 then
-
- put "infinite solution"
- goto program_exit
-
- end if
-
- exchange( row, col )
- update_array( row, col )
-
- end loop
-
- show_result
-
- program_exit:
-
- end program
-
- %
- % "initialize" sets up the simplex array for the constrained optimization
- % problem:
- %
- % maximize: y = 2*x1 - 3*x2 + 3*x3
- % subject to: 2 >= -2*x1 + 3*x2
- % 5 >= 1*x1 - 2*x2 - 4*x3
- % 6 >= 2*x1 + 1*x2 + 1*x3
- % 10 >= 1*x1 + 1*x2 + 5*x3
- %
- % and all variables are non-negative
- %
- procedure initialize
-
- % variable indices
- non_basis[0] := 0
- non_basis[1] := 1
- non_basis[2] := 2
- non_basis[3] := 3
-
- % slack variable indices
- basis[0] := 0
- basis[1] := N + 1
- basis[2] := N + 2
- basis[3] := N + 3
- basis[4] := N + 4
-
- % maximize y = 2*x1 - 3*x2 + 3*x3
- Sx[0,0] := 0
- Sx[0,1] := -2
- Sx[0,2] := 3
- Sx[0,3] := -3
-
- % 2 >= -2*x1 + 3*x2
- Sx[1,0] := 2
- Sx[1,1] := -2
- Sx[1,2] := 3
- Sx[1,3] := 0
-
- % 5 >= 1*x1 - 2*x2 - 4*x3
- Sx[2,0] := 5
- Sx[2,1] := 1
- Sx[2,2] := -2
- Sx[2,3] := -4
-
- % 6 >= 2*x1 + 1*x2 + 1*x3
- Sx[3,0] := 6
- Sx[3,1] := 2
- Sx[3,2] := 1
- Sx[3,3] := 1
-
- % 10 >= 1*x1 + 1*x2 + 5*x3
- Sx[4,0] := 10
- Sx[4,1] := 1
- Sx[4,2] := 1
- Sx[4,3] := 5
-
- end procedure
-
- %
- % "pivot_col" finds the column location for pivoting
- %
- function pivot_col : int
-
- var j, col : int
- var t : real
-
- col := 0
- t := 0
-
- for j := 1...N do
-
- continue when Sx[0,j] >= 0.0
-
- if t > Sx[0,j] then
-
- t := Sx[0,j]
- col := j
-
- end if
-
- end for
-
- return col
-
- end function
-
- %
- % "pivot_row" finds the row location for pivoting
- %
- function pivot_row( col : int ) : int
-
- var i, row : int
- var r, t : real
-
- row := 0
- t := 1.0e+99
-
- for i := 1...M do
-
- continue when Sx[i,col] = 0.0
-
- r := Sx[i,0] / Sx[i,col]
-
- if r > 0 and r < t then
-
- t := r
- row := i
-
- end if
-
- end for
-
- return row
-
- end function
-
- %
- % "exchange" switches the switches basis and non-basis variable indices
- % corresponding to the pivot row and column
- %
- procedure exchange( row, col : int )
-
- var t : int
-
- t := basis[row]
- basis[row] := non_basis[col]
- non_basis[col] := t
-
- end procedure
-
- %
- % "update_array" calculates the next set of values of the simplex array
- % for a given pivot row and column
- %
- procedure update_array( row, col : int )
-
- var pv, rv, cv : real
- var i, j : int
-
- % update pivot value
- pv := Sx[row,col]
- Tx[row,col] := 1 / pv
-
- % update pivot row element values
- for j := 0...N do
-
- continue when j = col
- Tx[row,j] := +Sx[row,j] / pv
-
- end for
-
- % update pivot column element values
- for i := 0...M do
-
- continue when i = row
- Tx[i,col] := -Sx[i,col] / pv
-
- end for
-
- % update the rest of the simplex array
- for i := 0...M do
-
- continue when i = row
-
- for j := 0...N do
-
- continue when j = col
- Tx[i,j] := Sx[i,j] - Sx[row,j] * Sx[i,col] / pv
-
- end for
-
- end for
-
- % copy temporary array into simplex array
- for i := 0...M do
-
- for j := 0...N do
-
- Sx[i,j] := Tx[i,j]
-
- end for
-
- end for
-
- end procedure
-
- %
- % "show_array" displays the simplex array
- %
- procedure show_array
-
- var i, j : int
-
- put " "...
-
- for j := 0...N do
-
- put non_basis[j]:12...
-
- end for
-
- put ""
-
- put " "...
-
- for j := 0...N do
-
- put "------------"...
-
- end for
-
- put ""
-
- for i := 0...M do
-
- put basis[i]:4, " |"...
-
- for j := 0...N do
-
- put Sx[i,j]:12...
-
- end for
-
- put ""
-
- end for
-
- put ""
-
- end procedure
-
- %
- % "show_result" shows the optimum value of the objective function
- % and the corresponding values of the variables
- %
- procedure show_result
-
- var i, j : int
- var found : boolean
-
- put "optimum y = ", Sx[0,0], " at:"
-
- for j := 1...N do
-
- put " x", j, " = "...
-
- found := false
-
- for i := 1...M do
-
- if basis[i] = j then
-
- put Sx[i,0]
- found := true
-
- end if
-
- end for
-
- if not found then
-
- put 0
-
- end if
-
- end for
-
- end procedure